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Abstract 

A partial differential equation model of a cantilevered beam with a tip mass at its free 
end is used to study damping in a composite. Four separate damping mechanisms consisting 
of air damping, strain rate damping, spatial hysteresis and time hysteresis are considered 
experimentally. Dynamic tests were performed to produce time histories. The time history 
data is then used along with an approximate model to form a sequence of least squares 
problems. The solution of the least squares problem yields the estimated damping coeffi- 
cients. The resulting experimentally determined analytical model is compared with the time 
histories via numerical simulation of the dynamic response. The procedure suggested here 
is compared with a standard modal damping ratio model commonly used in experimental 
modal analysis. 


1 This research was supported by the National Aeronautics and Space Administration under NASA Con- 
tract Nos. NAS1-18107 and NAS1-18605 while the first author was in residence at the Institute for Computer 
Applications in Science and Engineering (ICASE), NASA Langley Research Center, Hampton, VA 23665. 
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I. Introduction 

In this paper a variety for damping mechanisms for a quasi-isotropic pultruded composite beam 
are examined. The approach taken here is a physical one. The beam is modeled by a partial 
differential equation describing the transverse vibration of a beam with tip mass. The damping 
mechanisms considered are all physically based as opposed to the usual modal models. In total, 
four possible damping mechanisms are considered, one external and three internal. They are: 
viscous damping (air damping); strain rate damping; spatial hysteresis; and time hysteresis. 

In addition, various combinations of these mechanisms are considered. These physical damping 
models are incorporated into the Euler-Bemoulli beam equation, with care taken to formulate 
boundary conditions that are compatible with the various damping models. The resulting partial 
differentia] equation (integro partial differential equation in the case of time hysteresis damping) is 
approximated using cubic splines. The time histories of the measured experimental responses are 
then used to form a least squares fit-to-data parameter estimation problem. The mathematical 
arguments underlying this procedure are complete and imply convergence of a sequence of 
parameter estimates obtained from finite dimensional models to a set of best fit coefficients of the 
partial differential equation model. The least squares estimates of the various different damping 
parameters are then used in the partial (integro partial) differential equation to numerically simulate 
the response of the system. This numerically generated time response of the estimated system is 
then compared with the actual experimental time histories. These comparisons allow several 
conclusions to be drawn regarding the physical damping mechanisms present in the composite 
beam. 

In particular it is shown that the spatial hysteresis model combined with a viscous air damping 
model results in the best quantitative agreement with the experimental time histories. The results 
also support the physically intuitive notion that air damping should play a more significant role in 
lower modes while internal damping plays a more significant role for higher modes. It is also 
shown explicitly that the proposed damping models listed above cannot be modeled with any 
degree of success or consistency by using standard modal damping ratios, as the traditional modal 
analysis approach completely masks the physics of damping mechanisms. 

II. Basic Beam Model 

The beam considered here is a pultruded quasi-isotropic composite beam constructed for use in the 
proposed space station (Wilson and Miserentino, 1986). As such, the configuration of interest is 
a cantilevered beam with a mass attached to the free end. The beam is constructed of a biaxial 
(0°/90°) fiberglass roving held in place with knitted polyester yam with an equal volume of fibre 
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in both orientations. An isophtalic polyester resin system was used as the matrix. This material 
provides an alternative to aluminum which is lower in cost, has higher specific strength but is 
dynamically similar. As is illustrated here, this material also has interesting damping properties - 
dissimilar to those of aluminum. 

The equation of motion for the flexural vibration of a beam is easily calculated from consideration 
of the equilibrium of forces acting on a differential segment of beam (see for instance Clough and 
Penzien, 1975). In this formulation, damping can easily be included by adding the appropriate 
force or moment to the equations of equilibrium. A partial differential equation model of the beam 
with general damping is of the form 

+ L\U t (x,t. ) + L 2 u(x,t ) + ~ £ u xx (x,t)^ =f(x,t) (1) 

for x g (0, /), t > 0, subject to the appropriate boundary conditions and the initial conditions 
(taken to be u = u, = 0 at t = 0). Here p is the linear mass density (mass per unit length) of the 
beam, A is the cross sectional area of the beam, EI(x ) is the spatial varying flexural stiffness of 
the beam, the subscript indicates partial differentiation with respect to the indicated variable and 
u(x,t ) is the beam displacement in the transverse direction at location x, time t. The function u{x,t) 
is assumed to be smooth enough so that all the appropriate derivatives exist. The term L\U t {x,t) + 
L 2 u(x,t) is the focus of attention in this paper. The nature of the operator Lj is determined by 
external damping mechanisms while the nature of the term L 2 u(x,t ) is determined by internal 
damping mechanisms. 

The boundary conditions of interest here are those for a beam clamped at the end x = 0 and with a 
free end at x = /. In addition a mass of mass mj and rotational inertia J, is attached at x-l. The 
fixed end requires that the displacement and the slope of the displacement both be zero. This 
yields: 

«(0,r) = 0 (2) 

^(0,0 = 0 . (3) 

The free end requires that the sum of the moments at x = / and the sum of the forces acting at x - 1 
must each be zero. For die case of a tip mass at the free end, these boundary conditions become 
EI(l)u xx (l,t)=-Ju ta (l,t ) (4) 

[EHx)Uxx( x > l )} x = m T u u (l,t) , x =1 (5) 

as long as only external damping is present. 

Equations (1) - (5) describe the transverse vibration of a beam satisfying the Bemoulli-Euler 
assumption that the bending wave length is several times larger than the cross sectional 
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dimensions of the beam, and that only lower frequency excitations are applied to the beam. It is 
also assumed that rotary inertia of the beam, shear displacement of the beam and axial 
displacements arc negligible. 

If the tip mass is not present, the boundary conditions of equations (4) and (5) change 
accordingly. In addition, the nature of the damping operator L 2 will effect the boundary 
conditions. For the case of L\ = L 2 = 0, the vibration analysis problem is very simple as is the 
inverse problem addressed here. The nature of the damping mechanisms drastically changes the 
nature of the solution to the vibration problem and hence controls the response of the beam. The 
following section discusses several possible choices for modeling the operator L 2 in equation (1) 
and hence the internal damping mechanisms in the beam. 

III. Damping Models 

As mentioned in the introduction, four models for the damping mechanisms are examined. Two 
of these are time independent proportional models lending themselves to modal expansions, the 
other two are nonproportional hysteretic models. Various combinations of these models are also 
considered. 


Viscous Air Damping The most straight forward method of modeling the damping of a beam (or 
other object) vibrating in air is to use a viscous model with damping force assumed proportional to 
velocity. In this case the operator Lj becomes 

Lj = yl 0 (6) 

where 7 0 is the identity operator and y is the viscous damping constant of proportionality. The 
physical basis of this approach is a simple model of air resistance. As the beam vibrates it must 
displace air causing the force yu t (x,t) to be applied to the beam. Mathematically, this form of 
damping is used because it is proportional and easily treated using the same methods of analysis 
used for undamped systems (see Inman, 1989, for instance). This form of damping is often 
called external damping. 


Kelvin-Voigt Damping Kelvin- Voigt damping, or strain rate damping as it is sometimes called, is 
damping of the form 


L 7 = c d I 


35 

dx*dt 


(7) 


where 7 is the moment of inertia and c d is the strain rate damping coefficient. This model also 
satisfies a proportional damping criteria and hence is mathematically convenient It is compatible 
with theoretical modal analysis and is also, along with viscous damping, widely used in finite 
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element modeling. This form of damping is often referred to as internal damping and represents 
energy dissipated by friction internal to the beam. 

Unlike viscous external damping, inclusion of this form of damping affects the free end boundary 
conditions because it is strain dependent. The strain rate dependance results in a damping moment 
M d of the form 

M D = c d I(x) — — (8) 

ax l dt 

which is included in the derivation of the equation of motion (Clough and Penzien, 1975) and 
hence must also be included in any boundary conditions (such as a free end condition) depending 
on the moment. 

The full equation of motion and boundary conditions for the case including linear viscous external 
damping and Kelvin-Voigt internal damping can thus be written 

d 2 

p u tt (x,t)+ — [ElUja (x,t) + c d Iu^tixj)] + yu t (x,t)= f(x,t ) , xe (0,/), t > 0 

dx l 

u(0,r) = u x (0,t) = 0 , t>0 (9) 

Elliot) + c d I , t > 0 

— [Elu^x.t) + c d IUju^xj)] = m T u tt (x,t), x = l , t> 0. 
dx 

Here, note that the tip mass as well as the damping moment are represented in the boundary 
conditions The total damping mechanism used in (9) is the analog to proportional damping (i.e., a 
linear combination of mass and stiffness operators). 

Time Hysteresis Hysteretic damping terms are most commonly associated with sinusoidal 

loadings. The generic idea of including a mechanism in the beam vibration constitutive equation 

indicating that stress is proportional to strain plus the past history of the strain can be 

accomplished by introducing an integral term of the form 
0 

jg{s)u xx {x,t-^s)ds 

-r 

with the history kernel g(s) defined by 
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where a and |3 are constants. Since the introduction of the heredity integral occurs in the stress 
strain relationship, the boundary conditions must also be modified. In this case the boundary 
value problem of interest becomes 


d 2 

pU,l{x,t) + a^2 


V 

EIu xx (x,t) - J g(s)u xx (x,t+s)ds 


=f(x,t) ,xe (0,1), / >0 


u(0,t) = u x (0,t) = 0 t > 0 


( 12 ) 


Elu^ilt) - | g(s)u xx (l,t+s)ds - Ju xll (l,t), t > 0 


_a_ 

dx 


EIu xx (x,t) - J g(,s)M xx (;r,r+j)d.y 
-r 


= m T u„(x,t), x = l,t> 0 


It is emphasized again that the inclusion of a damping moment in the equation of motion also 
affects the boundary condition. 

Spatial Hysteresis Another type of damping (proposed by Russell, 1990) is based on interpreting 
the energy lost in the transverse vibration of a beam as resulting from differential rates of rotation 
of neighboring beam sections causing internal friction. This can be modeled by the damping 
expression 

d_ 

dx |.o’ 


r i 

jh(xX){u xt (t,x)-u xl (t^))d^ 


(13) 


where the kernel h(x£) may be defined, for example, by 


h(xft = 


.e -(x-tf/lb 2 ' 




(14) 


Under these circumstances the model for the beam vibration becomes (including viscous external 
damping) 

d 2 

P K tt (x,0 + r-r [Elu xx {x,t)] + y u,(x,t) 


d_ 

dx 


dx 2 

l 


Lo 


]h(x£)[u xt (x,t) - u x fy)}d^ 


=f(x,t), xe (0,/) , t > 0 


(15) 


u(0,t) = u x (0,t) = 0, r > 0 
EIu^ ( l,t ) = - Ju x[t ( l,t ) t > 0 
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= m T u tt (x,t ) , t> 0 , x = / 


/ 

\h(xX)[u xt (x,t) - u xl (^,t))dx 
0 

where again the internal damping mechanism is reflected in the boundary conditions. 

In total, the models described by systems (9), (12) and (15) represent four possible mechanisms 
of damping taken in various combinations. The approach taken here is to attempt to fit each of the 
combinations of damping models listed above to experimentally measured data. By examining 
each model’s numerical solution in comparison with measured data, a best model is chosen from 
these as being most representative of the cantilevered quasi-isotropic beam. As is discussed in 
Section V, these models all admit reasonable mathematical treatment. 

IV. Problem Statement 

The various damping coefficients introduced in the preceding discussion cannot be measured by 
static experiments. Thus, the damping constants y, c^, a, P, a, and b must all be estimated based 
on measurements taken from dynamic experiments. The procedure suggested here is to estimate 
various groups of damping parameters such as indicated in the three models of (9), (12) and (15). 
Once these coefficients are estimated they can be used in the model to produce a numerical 
simulation of the response of the structure under consideration subjected to identical experimental 
inputs. The analytical time response (with the estimated coefficients) is then compared to the 
experimentally measured time response. The model that best agrees with (predicts) the 
experimental response is then considered to be a valid and desirable physical model. 

In particular, several vectors of parameters q are defined, one for each model of interest. For the 
three cases discussed here they are: 

qi = [El, c d I, y] (16) 

which delineates the first damping model as defined by system (9). Here c 4 is the internal strain 
rate damping coefficient and y is the linear air damping coefficient The second model considered, 
as defined by system (12), is characterized by the parameter vector 

q 2 = [El, a, P] (17) 

where a and P characterize the time hystoretic damping term. The last model considered contains 
a combination of linear air damping, defined by the coefficient y, and spatial hysteresis defined by 

the constants a and b. The parameter vector for the third system defined by (15) is 

q 3 = [El, y, a ,b] . (18) 


d 

— [ElUjuixj)] - 
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Other combinations of the four damping mechanisms were considered but were dismissed as 
discussed in the later section on results. Even though the techniques (Spline Inverse Procedures) 
described below can readily be used to treat spatially varying coefficients El and c d I, the 
experiments for the efforts presented in this paper were performed on uniform beams. Hence, 
consideration in the this paper will hence forth be restricted to constant El and c d I. 

Note that in each case the parameter vector contains the flexural stiffness constant El. For most 
common materials El is tested, tabulated and well known. However in this case the material is a 
prototype composite with unknown material properties. Thus El is also estimated. Because of 
the relative size of the air damping coefficient c d , the term c d I is estimated. 

V. Mathematical Foundation of the Estimation Problem 

Two approaches to solving the problem of determining the coefficients in the vectors qi, q 2 and 
q 3 are formulated here. The first approach involves application of experimental modal analysis 
methods (see Ewins, 1988 or Inman, 1989) to a theoretical modal analysis of equation (1). This 
procedure, suggested by Clough and Penzien (1975), can only be applied to the problem of 
estimating qj because modal analysis is not applicable to the hysteresis terms in q 2 and q 3 . This 
modal approach is presented for comparison and because it represents a standard approach for 
measuring damping. However, modal approaches cannot be be used to solve the inverse problem 
for general constitutive elements. The inverse procedures suggested as the second approach here 
do not have this limitation and can be applied to systems with spatially dependent physical 
parameters (El, etc.) as well as exotic damping mechanisms. 

The second approach taken here is a nonmodal procedure applicable to all three estimation 
problems, and forms our proposed method. This method is based on a careful consideration of 
the distributed parameter nature of the test article. It consists of forming a sequence of finite 
dimensional approximations (Galerkin type cubic B-splines) to equation (1) with an associated 
least squares fit-to-data (see Banks and Kunisch, 1989, for a general discussion of these ideas). 
For each of the damping models presented in this paper, a corresponding sequence of approximate 
estimates qf can be shown to converge to a best-fit parameter qj for the original distributed 

parameter system (1) (or specifically for (9), (12), or (15)). Detailed convergence arguments are 
given in Banks et al (1983, 1986, 1988, 1989), Banks and Ito (1988). The computational 
algorithms proposed here are based on these considerations of the distributed parameter nature of 
the estimation problem and is called the Spline-based Inverse Procedure (SIP). 
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Modal Analysis (EMA) The typical experimental approach to measuring the damping in a 
structure is to use experimental modal analysis (EMA) to determine modal damping ratios and 
natural frequencies. These quantities can then be used to determine the physical parameters 
contained in the vector q j . In the case discussed here, the tip mass is removed (for ease in 
exposition) so that m T =7 = 0, and the unit interval is used (i.e., / = 1). The damping operators 
L] and L % become 


L 1 =“7 
pA 


, 2 =^ 4 - 

pA dx^dt 


(19) 


which commute with the stiffness operator so that a modal representation is possible (Caughey 
and O'Kelley, 1965). According to this theory (see for instance Inman, 1989) the solution of (9) 
can be written as a series of products of two functions, u m (x,t) = «m( 0 4> m (x), which satisfy (9) 
for each m and whose sequence of partial sums converges to the unique solution of (9). Here the 
normalized functions <\> m (x ) are the eigenfunctions (mode shapes) of the stiffness operator 

a 2 f ei a 2 ■) 


a ^ 2 


pA dx 2 


( 20 ) 


subject to the boundary conditions of system (9). These eigenfunctions satisfy (for EI constant) 


where = to^(p A/Ef) 2 and the <t> m satisfy the orthogonality condition 


( 21 ) 




(x)<t> m (x)= 5 nm . 


( 22 ) 


Here 5 ^ is the Kronecker delta and co^ are the undamped natural frequencies of the system. 


Substituting u m (x,t) into equation (1), multiplying by <j) m (jt) and integrating with respect to x over 
the interval (0,1), one finds that each a m (t) must satisfy 

1 

( y C A 1 n 4 V . . El 


<hn(*) + 


pA 


pA m 


ajt) + ~ Kt = Jf(xM m (x)dx, t > 0 (23) 


for all m = 1,2,3... . Equation (23) has a direct relationship to the frequency domain measured 
modal data available from EMA. The experimental modal analysis procedure assumes that the 
structure consist of some finite number of single degree of freedom oscillators of the form: 


'<(*) + zinAnflmb) + =/m 


(24) 
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where and & m are measured damping ratios and natural frequencies respectively. Comparing 
coefficients of a m and a m in equations (23) and (24) , one obtains 


and 


for each m. 


p 4 — = & 2 
Pm p A 


2 tA, 

p A pA m 


(25) 


(26) 


As outlined in Cudney and Inman (1989), the elastic modulus E may be estimated from equation 
(25) by a linear least squares fit by which one obtains 

K 47t 2 pA 


e 4 i 


m 


=i /p 


?■ 


(27) 


m 


where J m = to m /27t hertz with ( m = 1,2,... K, a given measured set. Equation (26) can 
be written down once for each pair in this measured set ( m = 1,2,... K, to produce a 


least squares determination of the damping parameters of the form 


Y 

-Cd- 


= B tz. 


(28) 


Here B t denotes the generalized inverse (least square) of the 2xK matrix 

"l |3 4 ,/" 

1 P \l 

a.-L . . 

P A . 

1 P^/_ 


(29) 


and z is the A^xl vector z = [2 2 ^^>2 2 t^Ad 7 ” of measured modal information. The 

entries in B are calculated from the analytical solution of the eigenvalue problem for the stiffness 
operator above with appropriate boundary conditions. 


Equation (28) (and/or a weighted version of it) can be used to estimate the distributed damping 
parameters for the problem involving qj. While intuitively obvious and straight forward, this 
modal based method requires that El must be constant and that P m must be known in closed form. 
In addition, this approach is not applicable to the noncommuting damping models involving q 2 
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and q 3 or to problems with spatially varying coefficients. This provides the motivation for the 
nonmodal approach developed next. 

Spline Inverse Procedure (SIP) An alternative to estimating q z - from measured modal data 
(frequency domain) is to formulate a parameter estimation problem based on measured time 
histories of the test structure's response. Let %(/,t,) denote the acceleration measurements at the 
tip of the beam (x = l) at various times q. The inverse problem of interest is then to find the vector 
of parameters q such that 

.. . £ . \ a .. . .2 (30) 

7(q) = £ \u u (l,t iy q) - %(/,*,-) I 

i=l 

is minimized where u(x, t, q) denotes the solution of equation (1) with the appropriate boundary 
and initial conditions corresponding to parameter values q. Here M is the number of tip 
acceleration measurements. 

This estimation problem cannot, of course, be solved analytically. However, an iterative 
optimization scheme coupled with an approximation method for the infinite dimensional system of 
equation (1) may be used. The procedure suggested here is outlined as follows. First equation 
(1) is approximated via Galerkin procedures using cubic spline elements (N is used to denote the 
approximation index) to yield an approximate finite dimensional version which is solved for u N (t). 
The approximate accelerations u (( are then used in the cost function of equation (30) to define the 

finite dimensional estimation problem of minimizing 

M 

J N ( q) = Y, - %('•';) I 2 (31) 

J==l 

The solution of this set of estimation problems yields a sequence of estimates, denoted {q N }, of 
best fit parameters. Under appropriate assumptions, this sequence is then shown to converge to 
q*, a vector of parameters for the fully distributed parameter model of equation (1) which 
minimizes /(q) of equation (30). The theoretical formulation of this approach is presented next. 
The required proofs are omitted but can be found in Banks and Ito (1988), Banks,et al (1983, 
1986, 1989), and Banks and Kunisch (1989). 

The SIP estimation algorithm is formulated in weak or variational form by multiplying equation 
(1) by \j/(x) and integrating over in the interval ( 0 , 1 ). This yields 

< u lt , y > + < Ljiq, \\r > + <L 2 U, \jr> = < f,\\r > (32) 

where the inner product < • , • > for the Hilbert space H = L 2 ( 0,1) is defined by 
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(33) 


l 

<<t) 1 ,<|> 2 > = |<t> ] (x)<t> 2 (*)<lx - 

Equation (32) must hold for all \\r in V, a Hilbert space continuously and densely imbedded in the 
Hilbert space H , containing the solution of (1) subject to the appropriate boundary conditions. 
(In this case V = {ye \j/(0) = \y'(0) = 0) where H 2 is the Sobolev space of functions 

possessing first and second derivatives in Z^O,/)). The terms in equation (32) are now identified 
with elements ii(t;q), u(r,q) and u(r,q) which evolve in time, and satisfy (in a generalized or weak 
sense - see Banks and Ito (1988) and Banks and Kunisch (1989)) the evolution equation 

u (r, q) + B(q) u(r,q) + A(q) u(r,q) =f(t), t> 0 (34) 

subject to the appropriate initial conditions. Here the explicit dependence on the parameter vector 
q is emphasized, while the solution of equation (34) is a function of time for each x (or, 
alternatively, is thought of as a function of x for each t, the function being an element of V for 
each t). The operators A and B can be appropriately defined using the corresponding terms in 
(32). 

The Galerkin approach employing cubic spline subspaces to solve (34) (or equivalently (32)) 
iteratively is explained next. Given a value of N and a vector q, an approximate solution to (34) 
in X N = span { B B N ) is sought of the form 

N N 

w%-q) = ^ "i')B V = Xw"(r; q)Z^ (35) 

J = 0 j=0 

where { B ^ ) is the set of cubic spline basis functions appropriately modified to be in the domain 

KJ ft 

of defnition of the operators in equation (34). More precisely, let Ar = {x/} j= q with = il/N for i 

= 0, 1,..., N, and let B j,j = - 1 ,..., N + 1 denote the standard C 2 (0,/) basis elements for the 

cubic B-spline subspaces of dimension N + 3 corresponding to the grid A N (see Prenter, 1975). 
Here C 2 (0,1) is the set of all continuous functions with continuous first and second derivatives on 
the interval (0,/). Then B ^ is given by 

B^ = Bj, for 2 <7 < A/ - 2 
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<=< 4 ?;. 




Note then that X N = Sq (A w ) = (4>e S 3 (A^): 4>(0) = ())(/) = 0} where .S 3 (A w ) = {<()e C 2 (0,/); <)) is a 
cubic polynomial on each interval [x,vx j+ i]}. 

The appoximate solutions to (34) are determined from requiring that for all functions zeX N 

<u N (t), z>+ < B(q)w N (t), z > + <A(q) u N (t), z > = </(r), z > (36) 

with appropriate projected initial condition. Choosing z = B ^ and using equation (35) one may 
write this in an equivalent matrix form as 

w* (0 + D n \v N (t) + = F^(/) (37) 

where is the N+lxl vector [h^, and where the N+lxN+1 "damping" and 

"stiffness" matrices D N and K N are defined by 

Z^ = <fif,B(q)^> 

, A(q )B"> 

Here the subscript ij denotes the if* 1 element of the matrix and the vector is defined as the Nx 1 
vector of elements F? = c f(t), B^ >. Each approximate identification problem now reduces to 

calculating q^ that minimizes (31) subject to the vector differential equation (37). These 
calculations result in the sequence {q^} which, as mentioned above, converges to a q minimizing 
(30) subject to systems (9), (12) or (15), as appropriate. 

It is important to note here that this approximation differs from a standard finite element methods 
in two fundamental ways. First, the damping mechanism produced by the matrix D N converges 
to a physical damping model. Typical finite element models treat damping in an ad hoc fashion. 
Secondly, the entire model converges to a strength of materials/continuum mechanics model 
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having more physical significance than the rather arbitrary node model produced by standard finite 
element approximations. 

VI. Test Procedures 

This section describes a series of experiments performed to examine the proposed damping 
mechanisms and validate the approach (SIP) proposed in this paper. Figure 1 is a schematic of 
the experimental configuration. The experiments were performed in the Mechanical Systems 
Laboratory at the State University of New York at Buffalo. The beam is set up in a cantilevered 
arrangement with a removable and adjustable tip mass. The geometric properties of the beam are 
listed in Table 1. 

Substantial effort was made to insure that the fixed end behaved as a true clamped boundary 
condition. The clamp consisted of a large iron mass held together with a very strong magnetic 
field. The beam was extended an equal length beyond the clamp. The motion of the clamp and the 
extended portion of the beam were monitored to insure that no energy was lost to motion of the 
clamp and that no motion was transmitted through the clamp to the extended beam. The 
theoretical values of for a clamp-free beam (i.e., the solution to the transcendental characteristic 

equation) are listed in Table 2 for reference. 

The beam was excited by hammer hits. The hammer used was a Kistler Instrument Company 
(KIC) model with a load cell at its tip which allows a measurement of the input force. The tip of 
the hammer has variable mass and stiffness to insure excitation of the beam in the frequency range 
of interest. The hammer's load cell signal was run through a signal conditioning device (KIC) to 
the input channel of a GENRAD 2515 data acquisition and analyzer system for processing. 

It should be noted that many experimentalists suggest using a shaker along with swept sine or 
random excitation so that time histories of the forcing functions are repeatable in subsequent tests. 
However, because of the flexible nature of the structure (the first mode natural frequency is 3.64 
Hz) it was observed that attaching shaker to the structure either at the root or near the tip of the 
beam altered the measured natural frequencies. Consequently impact tests were used with the 
exception that the first mode was confirmed by analyzing the free response of the beam to an 
initial deflection in the first mode shape. 

The beam response was measured by using KIC piezobeam accelerometers and associated signal 
conditioning as well as a laser vibrometer (DISA), which provides a direct velocity measurement 
and a check on the accelerometer signals. The output responses are fed to the remaining channels 
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of the GENRAD analyzer. Once the response signals are in the analyzer they are digitized and 
stored for various programs written to perform the above analyses. 

In the modal analysis approach, the digital records of the input and output signals were 
manipulated into transfer functions and analyzed in the frequency doman using the Nyquist circle 
fit method (see Ewins 1984, for instance). The computed modal data and 0) m was averaged at 
least 15 times for each mode. The EMA program used was a standard commercially available 
EMA package (SDRC) and the results were checked using several other standard frequency 
domain algorithms. 

One popular alternative to frequency domain EMA is to use time domain realization methods. 
This approach uses the digitized time histories to construct a finite dimensional state matrix for the 
test structure. The numerical solution of the eigenvalue problem for this state matrix then yields 
the desired modal constants t, m and & m . One algorithm for this is the Eigensystem Realization 
Algorithm (ERA) developed at NASA by Juang and Pappa (1985). Data from the analyzer was 
transported by direct link to a VAX 1 1785 containing ERA. The results obtained by this approach 
agreed with those of the Nyquist approach lending credibility to the measured modal data. 

The experimentally determined modal data was then down loaded to a commercial matrix 
computation package (Moler et al 1987) where various least squares algorithms were used to solve 
equation (28) for the desired distributed parameters y and Cj. 

The time domain test data was also sent via BITNET to Brown's IBM 3090. The above described 
SIP methodology was applied to this data and qi, q 2 and qj estimated. Several tests were run in 
various configurations to investigate the effects of various tip mass arrangements as well as to 
examine various combinations of damping mechanisms. A large number of tests were performed 
over a two year period; in the following section, the findings of these tests are summarized. 

VII. Results 

The results of estimating the various damping parameters using EMA are discussed first as it is 
limited to the problem defined by qj. In the EMA approach, the correctness of a given estimate is 
judged by the ability of the solution to produce frequency independent parameters y and cj. Next, 
the results of the SIP approach are used to examine the damping mechanism. In this case the 
success of a given estimate is judged on the ability of the estimate to numerically simulate the 
experimental time history of the structure’s response. 
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Experimental Modal Analysis The results of performing 15 modal tests, as outlined in Section VI, 
are listed in Table 3, for the system without a tip mass. Note the large damping ratios exhibited 
by this material when compared with calculated values of aluminum or steel. This data was first 
used in equation (27) to determine the values of the modulus (£) calculated to be 

£ = 2.68 x 10 10 N/m 2 (38) 

with a variance of 0.6 N/m 2 . 


The excellent fit provided by the Bemoulli-Euler beam stiffness to the measured frequencies 
indicates that this is a suitable stiffness model for this particular composite with 0°/90° 
orientations. Some researchers have suggested that a Timoshenko model might be more 
appropriate for composites. However the inclusion of rotary inertia and shear effects did not 
provide a convincing fit to the data obtained here. 

It has been shown by Cudney and Inman (1988) that attempting to use just air damping, y, or just 
strain rate damping, c d , alone fails to match the measured modal data. In each case the attempt to 
fit a single damping parameter is measured by the ability of the estimated values of y and c d , to 
reproduce the measured damping ratios t, m . This situation is discussed later in the context of the 
SIP estimates. The significance of this result is that a single modal damping ratio cannot logically 
be used to model the damping mechanism of the composite beam. 


Next the generalized inverse of the data matrix B defined by equation (29) and the theoretical 
values of the eigenvalues given in Table 2 are used to calculate the desired damping coefficients y 
and c d from equation (28) and the vector z. The vector z contains the experimentally determined 
modal data of Table 3. The results of y and c d are 


y = 1 .7561 N-sec/m 2 , c d = 2.05 x 10 5 N-sec/m 2 (39) 

for 9 modes of data. 


The effect of natural frequency on the measured modal damping ratio is seen by substituting 
the analytical expression for into equation (26). This yields that 




2co 


m 


c d 

2E 2 / 


'rrr 


(40) 


This indicates clearly that the effect of air damping (8) decreases with increasing mode number 
(co^ — > °°). Thus for higher modes the strain rate damping makes a more significant contribution 
to the measured damping ratio. This agrees with the physically intuitive notion that the low 


16 



frequency modes are pushing more air than the higher frequency lower amplitude modes. In fact, 
for a free-free configuration it is claimed by Vinson (1989) that the effect of air damping can be 
subtracted based on Blevins' equation (Blevins, 1977) which considers flow effects. 

To check the validity of this model, the frequency dependence of the coefficients y and was 
examined by recalculating them using a different number of modes. Successive least squares was 
performed using first 2 modes, then 3 modes, etc., up to the total of 9 modes. The result is 
illustrated in Table 4. Table 4 indicates that the estimates of y and Crf depend somewhat on the 

frequency range of interest. This is inconsistent with the physical models put forth in Section HI 
for estimating q j. Furthermore this approach is not applicable to the problem of estimating q 2 and 
q 3 . Hence, one must conclude that the modal approach is not satisfactory in attempting to model, 
in composite beams, any of the damping mechanisms proposed in this paper. 

Spline Inverse Procedure The use of the SIP provides a nonmodal approach appropriate for each 
of the problems of estimating qj, q 2 and q 3 . The problem of estimating qj is solved first for 
comparison with the modal approach. This problem was solved again using a slightly more 
complicated cantilevered beam with a tip mass. The stiffness parameter (elastic modulus) E was 
estimated to be 2.71 x 10 10 N/m 2 , in good agreement with the modal estimation results above. 
Estimates of air damping alone or strain rate alone as a damping model proved to be inadequate in 
reproducing time histories matching those of the experimental data, indicating a poor model. This 
is again in agreement with the results obtained by using the modal approach. 

The time domain approach provided by SIP allows a convenient comparison between the 
measured time response and the analytical time response generated by the model of equation (9) 
with the experimentally determined parameter vector qj. The difference between the numerical 
solution for the time history of the acceleration u tt (Xj,t) for the analytical model with the estimated 
parameter qj and the experimentally measured accelerations define the residual which is generally 
small (Banks et al., 1987). The analytical time response is plotted along with the measured time 
response versus time in Figure 2. While the agreement is fair, the residual is larger for longer 
time intervals, warranting further modeling. 

Next the temporal hysteresis model involving q 2 was considered as a possible candidate for 
modeling the damping in the composite. In this case, the estimation procedure produces (i.e., 
consistent with our previous methods for estimating qj) a good value for E but drives the air 
damping coefficient to zero. The residual, however, is better than that for the model with qj. 
Figure 3 illustrates a plot of the measured acceleration versus time as well as the acceleration 
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predicted by the model with the estimate q 2 - The difference between the measured and predicted 
value over the time interval of interest is almost negligible. Because this model drives the air 
damping coefficient to zero (violating physical intuition), a third model (q 3 ) was considered. 

The last model considered is based on a concept of spatial hysteresis as defined by the estimation 
problem for qj. Again the resulting estimate of the elastic modulus E is consistent with those 
estimated previously. The values estimated for the spatial hysteresis parameters (a = 1.040394, 
b = 0.064362) and an air damping coefficient (y = .090189) produce an excellent match between 
predicted and measured response as indicated in Figure 4 (Banks, et al. 1988). However, the 
external damping coefficient y differs from that estimated using the parameter vector qj (y =.0315) 
emphasizing the fact that air damping should not be estimated independently whenever internal 
damping mechanisms are present. 

VIII. Conclusion and Discussion 

Three different models of damping have been presented to account for the experimentally 
observed dissipation in a pultruted composite beam. A spline based inverse procedure (SIP) 
which relies on the distributed parameter nature of the damping mass and stiffness parameters was 
proposed and used to estimate the form of each damping mechanism. External air damping, strain 
rate damping, spatial hysteresis and time hysteresis models were considered. The spline based 
method was also compared to a standard experimental modal analysis (EMA) approach. The 
EMA approach is not applicable to the various hysteresis models, nor is it applicable to systems 
with spatially varying parameters in general. Both the SIP and EMA approaches yield consistent 
values for the elastic modulus ( E) for all three estimation models. This is consistent with the fact 
that frequencies are much more robust to estimates than are damping quantities. 

Both hysteresis models produce better results than the strain rate damping model. However, the 
spatial hysteresis model allows for the air damping term while time hysteresis does not. Since air 
damping is obviously present, the time hysteresis result is less satisfying. A comparison of the 
hysteresis models is given in Banks, et al (1988). As indicated in that presentation further 
analysis and modeling is required before a conclusive decision can be made about a best model. It 
is clear from the results presented here that hysteretic damping is able to reproduce experimental 
time responses with more accuracy than the standard Kelvin-Voigt model. It is also clear that the 
standard method of measuring damping, EMA, does not provide an accurate method for 
investigating the energy dissipation in the composite beam tested here. 
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In summary, a new method of determining damping mechanisms in a distributed parameter model 
has been proposed and applied to a beam. The method refereed to as SIP, also yields estimates of 
the stiffness parameters. This method has been compared to the standard method of determining 
damping in structures using modal methods. When compared on the same experimental test data, 
the SIP approach produces more consistent estimates of the Kelvin-Voigt damping parameter then 
those obtained by using modal methods. In addition, the proposed procedure is applicable to 
hysteretic damping models and to systems with spatially varying parameters that cannot be treated 
by modal methods. 
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Table 1. Beam parameters. 


Length (meters) (1) 

1.0 

Moment of Inertia (meter 4 ) (/) 

1.64 x lO' 9 

Density (kilograms/meter 3 ) (p) 

1710 

Area, cross section (meter 3 ) (A) 

0.597 x 10- 3 


Table 2. Theoretical eigenvalues of a clamped-free beam. 


Mode Number 

Eigenvalue (J3j) 

1 

1.875 

2 

4.694 

3 

7.855 

4 

10.996 

5 

14.137 

6 

17.279 

7 

20.420 

8 

23.562 

9 

26.704 


Table 3. Experimentally measured modal data. 


Mode 

Theoretical 
Freq. (Hz.) 

Experimental 
Freq. (Hz.) 

Damping 
Ratio (%) 

Std. 

Dev. 

2 

23.096 

22.8 

.218 

.015 

3 

64.675 

65.3 

.227 

.016 

4 

126.740 

127 

.154 

.003 

5 

209.488 

212 

.228 

.023 

6 

312.955 

314 

.120 

.008 

7 

437.075 

435 

.131 

.021 

8 

581.475 

580 

.155 

.010 

9 

747.475 

733 

.202 

.015 
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Table 4. Estimates of damping based on using a successive 
number of nodes (all values in N sec/m2). 


Modes 

Viscous Only 
Ci 

Strain Only 
Q.xlO 6 

Viscous and Strain 
C] C^xlO^ 

1-2 

.3619 


.0724 


1-3 

.8755 

.3127 

.2014 

.2699 

1-4 

1.2829 

.1179 

.6053 

.0873 

1-5 

2.2693 

.0978 

.6157 

.0856 

1-6 

2.6962 

.0451 

1.3901 

.0323 

1-7 

3.3565 

.0304 

1.6867 

.0221 

1-8 

4.3778 

.0251 

1.8039 

.0199 

1-9 

6.0027 

.0236 

1.7561 

.0205 



Figure 1. Schematic of test configuration. 
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